﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace ExportToDWGTool
{
    class ConvertPlaneCoord
    {
        /**	 * 平面坐标系	 */
        private static double a = 6378137; // 赤道圆的平均半径	
        private static double f = 1 / 298.257223563; // WGS84椭球扁率		
                                         	
        public static void turn(double lLongtitude, double bLatitude,double h,ref double x,ref double y,ref double z)
        {         
            double b = a * (1 - f); //短半轴

            double e = Math.Sqrt(a * a - b * b) / a;  //偏心率

            b = 6356755.00;
            e = 0.016710219;

            //曲率半径
            double N = a / Math.Sqrt(1 - e * e * Math.Sin(bLatitude * Math.PI / 180) * Math.Sin(bLatitude * Math.PI / 180));

            double WGS84_X = (N + h) * Math.Cos(bLatitude * Math.PI / 180) * Math.Cos(lLongtitude * Math.PI / 180);

            double WGS84_Y = (N + h) * Math.Cos(bLatitude * Math.PI / 180) * Math.Sin(lLongtitude * Math.PI / 180);

            double WGS84_Z = (N * (1 - (e * e))) * Math.Sin(bLatitude * Math.PI / 180);

            x = WGS84_X;
            y = WGS84_Y;
            z = WGS84_Z;

        }


        ///<summary>
        /// GPS经纬度换算成x,y坐标
        ///</summary>
        ///<param name="l">经度</param>
        ///<param name="B">纬度</param>
        ///<param name="xc">X坐标</param>
        ///<param name="yc">Y坐标</param>
        public static void MCT84Bl2xy(double l, double B, out double xc, out double yc)
        {
            try
            {
                l = l * Math.PI / 180;
                B = B * Math.PI / 180;

                double B0 = 30 * Math.PI / 180;

                double N = 0, e = 0, a = 0, b = 0, e2 = 0, K = 0;
                a = 6378137;
                b = 6356752.3142;
                e = Math.Sqrt(1 - (b / a) * (b / a));
                e2 = Math.Sqrt((a / b) * (a / b) - 1);
                double CosB0 = Math.Cos(B0);
                N = (a * a / b) / Math.Sqrt(1 + e2 * e2 * CosB0 * CosB0);
                K = N * CosB0;

                double Pi = Math.PI;
                double SinB = Math.Sin(B);

                double tan = Math.Tan(Pi / 4 + B / 2);
                double E2 = Math.Pow((1 - e * SinB) / (1 + e * SinB), e / 2);
                double xx = tan * E2;

                xc = K * Math.Log(xx);
                yc = K * l;
                return;
            }
            catch (Exception ErrInfo)
            {
            }
            xc = -1;
            yc = -1;
        }
    }
}
